Method of determining a map of height of liquid hydrocarbon in a reservoir

ABSTRACT

This invention relates to a method for determining a map of expectation and/or of variance of height of liquid hydrocarbon in a geological model. The method allows analytical resolution of these maps while taking account of the uncertainties in the variables allowing this calculation such as the porosity or oil saturation of the rock and also the uncertainty in the presence of certain types of facies given by the apportionment cubes of architectural elements p AE (c) and proportion cubes for each facies, these proportions then having a triangular distribution defined by the following three values p A,AE,min (c), p A,AE,max (c) and p A,AE,mode (c). In particular, the method comprises the calculation of the sum of the plurality of architectural elements of p AE (c)·⅓(p A,AE,min (c)+p A,AE,max (c)+p A,AE,mode (c)) and the determining of a value of expectation of height of liquid hydrocarbon for said column as a function of the sum determined.

RELATED APPLICATIONS

The present application is a National Phase entry of PCT Application No.PCT/FR2015/051309 filed May 19, 2015, which claims priority from FRPatent Application No. 14 54477, filed May 19, 2014, said applicationsbeing hereby incorporated by reference herein in their entirety.

FIELD OF THE INVENTION

This invention relates to the field of determining liquid hydrocarbonreserves in a reservoir and, in particular, the field of maps of heightof liquid hydrocarbon in a reservoir.

BACKGROUND OF THE INVENTION

In order to easily determine the hydrocarbon reserves (gas, petrol oroils) of reservoirs in operation or reservoirs to be operated,geologists most often manipulate and study two-dimensional maps thatrepresent the “HuPhiSo” (referred to as “HuPhiSo” maps). Indeed, thesemaps are summarised and quickly provide the geologist with informationas to the economic value of the reservoir (3D).

For a three-dimensional rock volume, the “HuPhiSo” map (in twodimensions) comprises a plurality of values, with each value being theproduct of several parameters:

HuPhiSo=h·NTG·φ·So

with the parameter h being the height of the column (according to theaxis {right arrow over (z)} if the map is according to the plane ({rightarrow over (x)},{right arrow over (y)})), the parameter NTG (“Net toGross”) being the fraction of the rock volume favourable to the presenceof hydrocarbons, with the parameter φ being the porosity of the rock andthe parameter So being the oil saturation of the porosities.

As such, if V is the rock volume considered, V·NTG·Phi So is the volumeof oil contained in the rock volume. A value “HuPhiSo” represents the“height of the column of oil” in a column under consideration. However,this value is indirectly connected to the volume of the oil and canrepresent the “surface density of the oil”.

For a model meshed with parallelepiped cells, it is possible tocalculate a “HuPhiSo” value for the final map by summing the value“HuPhiSo” of each cell of a column of this model.

Of course, the determining of the “HuPhiSo” values can be simple if theparameters h, NTG, Phi and So are known for certain for any point/cellof the model. However, in practice, these parameters are perfectly knownonly in a few isolated points: the points of the drilling wells. Ofcourse, it is possible to simulate these parameters for the other cellsthanks to the crossing of various data analyses (ex. seismic) and theuse of stochastic simulation algorithms based on geostatistical methodsin order to “supplement” the model. To this is added the naturaldistribution of the petrol-physical properties (acoustic impedance,effective porosity, porosity of the rocks in the NTG, etc.). Eachphysical magnitude can therefore be considered as a random variable,possibly non-stationary (ex. the porosity decreases most often with thedepth).

There are a large number of uncertainties that weigh on the variousparameters on which the HuPhiSo map depends, in such a way that thelatter can vary substantially from one possible map to another.

For making decisions, reservoir engineers particularly appreciate themaps of expectation of the HuPhiSo and maps of expectation (or ofvariance) around the expectation of the HuPhiSo.

Usually, these maps are obtained by simulating a large number of HuPhiSomaps. For each cell of the map, an average and a variance are thencalculated according to the corresponding cells in the simulated HuPhiSomaps.

However, each simulation can be costly in terms of resources and incalculation time and it is not rare for the determining of maps ofexpectation and standard deviation to take several hundred hours.

There is therefore a need for determining an average HuPhiSo map(expectation) and/or of the interval of uncertainties corresponding tothis map (variance or standard deviation) more quickly, without carryingout a large number of simulations of HuPhiSo maps.

SUMMARY OF THE INVENTION

To this effect, this invention proposes an analytical method that makesit possible, after certain simplifications, to directly determine theexpectation of the height of liquid hydrocarbon in a geological model.

This invention then aims for a method of determining a map ofexpectation of height of liquid hydrocarbon in a geological model, withthe geological model comprising meshes able to be associated with afacies among a plurality of facies and with an architectural elementamong a plurality of architectural elements.

The method comprises, for at least one column of the geological model:

/a/ receiving an apportionment cube for each architectural element AEamong the plurality of architectural elements, with the apportionmentcube describing for each mesh c of the model the probability p_(AE)(c)that said mesh belongs to said architectural element,

/b/ receiving a proportion cube for each facies among the plurality offacies and for each architectural element among the plurality ofarchitectural elements, with the proportion cube associated with afacies A and with an architectural element AE describing for each mesh cof the model the proportion of said facies in said mesh if said meshbelongs to said architectural element AE, with said proportion being arandom variable of triangular distribution defined by a minimum valuep_(A,AE,min)(c), a maximum value p_(A,AE,max)(c), and a modep_(A,AE,mode)(c),

/c/ determining, for each mesh c of said column and for each facies A ofthe plurality of facies, the sum over the plurality of architecturalelements of p_(AE)(c)·⅓(p_(A,AE,min)(c)+pA,AE,maxc+pA,AE,modec.

/d/ determining of the value of expectation of height of liquidhydrocarbon for said column analytically as a function of the sumdetermined in the step /c/ for each mesh c of said column and for eachfacies A.

Normally, the expectation of the value HuPhiSo is very difficult tocalculate directly Esp(HuPhiSo)=Esp(h·NTG·φ·So) because all of thevariables NTG, Phi So depend in particular on the facies and as such onthe location of the facies in the model.

The mode of a triangular distribution can also be called base.

In this invention, Esp(X) designates the expectation (average value) ofa random variable X, Var(X) designates its variance and Cov(X,Y)designates the covariance of two random variables X, Y.

By making the hypothesis of the independence of the parameters NTG, Phiand So, it is possible to greatly simplify the calculation of thisexpectation. Experimentally this hypothesis was verified and thissimplification, which is not the simple application of mathematicalprinciples, makes it possible to mutualise the calculation of the sum(over the plurality of architectural elements) of p_(AE)(c)·⅓(p_(A,AE,min)(c)+p_(A,AE,max)(c)+p_(A,AE,mode)(c)) and possible toseparately calculate the expectations of Esp(φ_(A)(c)), Esp(NTG_(A)(c)),Esp(So_(A)(c)).

This simplification allows for better effectiveness during thecalculations to be made of the expectation and as such to substantiallyimprove the required calculation time.

Furthermore, with each mesh c of the model able to be associated with aproportion of porosity for each one of the facies of said plurality offacies, with said proportion of porosity able to be a random variable oftriangular distribution defined by a minimum value φ_(A,min)(c), amaximum value φ_(A,max)(c), and a mode φ_(A,mode)(c), the method canfurther comprise:

/c-1/ determining, for each mesh c of said column and for each facies Aamong the plurality of facies, the value⅓(φ_(A,min)(c)+φ_(A,max)(c)+φ_(A,mode)(c)),

The determination of the step /d/ can then be according to the valuedetermined in the step /c-1/ for each mesh c of said column and for eachfacies A among the plurality of facies.

As such, the value Esp(φ_(A)(c)) can then be calculated easily,independently of the other parameters, if the variable φ_(A)(c)comprises an uncertainty.

In a possible embodiment, with each mesh c of the model able to beassociated with a proportion of rock volume favourable to the presenceof hydrocarbons for each one of the facies of said plurality of facies,with said proportion of rock volume favourable to the presence ofhydrocarbons able to be a random variable of triangular distributiondefined by a minimum value NTG_(A,min)(c), a maximum valueNTG_(A,max)(c), and a mode NTG_(A,mode)(c), the method can furthercomprise:

/c-2/ determining, for each mesh c of said column and for each facies Aamong the plurality of facies, the value⅓(NTG_(A,min)(c)+NTG_(A,max)(c)+NTG_(A,mode)(c)),

The determination of the step /d/ can then be according to the valuedetermined in the step /c-2/ for each mesh c of said column and for eachfacies A among the plurality of facies.

As such, the value Esp(NTG_(A)(c)) can be calculated easily,independently of the other parameters, if the variable NTG_(A)(c)comprises an uncertainty.

Advantageously, with each mesh c of the model able to be associated witha liquid hydrocarbon saturation, with said liquid hydrocarbon saturationable to be characterised by a liquid hydrocarbon-water interfacedimension z_(w) that has a probability p(z_(w)) of distribution betweena maximum value z_(w,min), and a minimum value z_(w,max), and by atriangular distribution of the hydrocarbon saturation for each dimensionh above this interface dimensions, with said distribution of hydrocarbonsaturation able to be defined by a minimum value CS_(HC,A,min)(h), amaximum value CS_(HC,A,max)(h), and a mode CS_(HC,A,mode)(h) the methodcan further comprise:

/c-3/ determining, for each mesh c of said column having a dimension zand for each facies A among the plurality of facies, the value

$\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{{cs}_{{HC},\max}\left( {z - z_{w}} \right)} + {{cs}_{{HC},{mode}}\left( {z - z_{w}} \right)} + {{cs}_{{HC},\min}\left( {z - z_{w}} \right)}}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}$

The determination of the step /d/ can be a function of the valuedetermined in the step /c-3/ for each mesh c of said column and for eachfacies A among the plurality of facies.

By notation, it is chosen that the notation z_(w,min) corresponds to thehighest dimension (with an axis {right arrow over (z)} upwards) incomparison with the dimension z_(w,max). As such, we havez_(w,min)≧z_(w,max) by taking this notation. Indeed, the dimension ofthe interface z_(w,min) (respectively z_(w,max)) corresponds to thedimension for which the quantity of hydrocarbon is the least substantial(respectively the most substantial).

All of the notations of dimension z_(X,min), and z_(X,max) respect thesame principle (with X able to be any character string).

In a particular embodiment of the invention, with said liquidhydrocarbon saturation able to be characterised by a dimension of theliquid hydrocarbon-gas interface z_(g) having a probability p(z_(g)), ofdistribution between a minimum value z_(g,min), and a maximum valuez_(g,max), the method can further comprise:

/c-4/ determining, for each mesh c of said column having a dimension z,the value

$\left\{ {1 - {\int_{z_{g,\max}}^{z}\begin{matrix}{{0\mspace{14mu} {if}\mspace{14mu} z} \geq z_{g,\min}} \\{{{p(z)}_{g}{dz}_{g}\mspace{14mu} {if}\mspace{14mu} z_{g,\max}} \leq z \leq z_{g,\min}} \\{{1\mspace{14mu} {if}\mspace{14mu} z} \leq z_{g,\max}}\end{matrix}}} \right.$

The determination of the step /d/ can be a function of the valuedetermined in the step /c-3/ for each mesh c of dimension z of saidcolumn and for each facies A among the plurality of facies, multipliedby the value determined in the step /c-4/ for each mesh c of dimension zof said column.

Furthermore, this invention can also relate to a method for determininga map of variance of height of liquid hydrocarbon in a geological model.

The method can comprise, for at least one column of the geologicalmodel:

/e/ determining of a value of expectation of height of liquidhydrocarbon for said column, as disclosed hereinabove;

/f/ determining of a value of variance of height of liquid hydrocarbonfor said column analytically according to the sum determined in the step/e/.

In an embodiment, the values of the proportion cubes for two differentcells of the model and for different architectural elements can beconsidered independently between them.

Furthermore, the method can further comprise:

/g/ receiving a correlogram ρ(A, Δz) according to a direction of saidcolumn for said facies A in a given architectural element;

/h/ determining the expectation of the product of the presence of thefacies A in said architectural element in a mesh c1 by the presence ofthe facies A in said architectural element in a mesh c2, with thedistance between c1 and c2 according to said direction of thecorrelogram being Δz, with the probability of the presence of the faciesA of said architectural element for the mesh c1 being p_(A,AE,c1), withthe probability of the presence of the facies A of said architecturalelement for the mesh c2 being p_(A,AE,c2), according to:

$\left( {1 - {\rho \left( {A,{\Delta z}} \right)}} \right) \cdot \left( {{{{Esp}\left( p_{A,{AE},{c1}} \right)} \cdot \left( {{{Esp}\left( p_{A,{AE},{c2}} \right)} + \sqrt{{{Var}\left( p_{A,{AE},{c1}} \right)} \cdot {{Var}\left( p_{A,{AE},{c2}} \right.}}} \right)} + {{\rho \left( {A,{\Delta z}} \right)} \cdot \left( {{\frac{1}{2}{{Esp}\left( p_{A,{AE},{c1}} \right)}} + {\frac{1}{2}{{Esp}\left( p_{A,{AE},{c2}} \right)}}} \right)}} \right.$

The determination of the step /f/ can then be according to thedetermination of the step /h/.

This invention also relates to a device intended for determining a mapof expectation of height of liquid hydrocarbon in a geological model,with the geological model comprising meshes able to be associated with afacies among a plurality of facies and with an architectural elementamong a plurality of architectural elements.

The device comprises, for at least one column of the geological model:

/a/ an interface for receiving an apportionment cube of architecturalelement for each architectural element AE among the plurality ofarchitectural elements, with the apportionment cube describing for eachmesh c of the model the probability p_(AE)(c) that said mesh belongs tosaid architectural element,

/b/ an interface for the receiving of a proportion cube for each faciesamong the plurality of facies and for each architectural element amongthe plurality of architectural elements, with the proportion cubeassociated with a facies A and with an architectural element AEdescribing for each mesh c of the model the proportion of said facies insaid mesh if said mesh belongs to said architectural element AE, withsaid proportion being a random variable of triangular distributiondefined by a minimum value p_(A,AE,min)(c), a maximum valuep_(A,AE,max)(c), and a mode p_(A,AE,mode)(c),

/c/ a circuit for the determining, for each mesh c of said column andfor each facies A of the plurality of facies, the sum over the pluralityof architectural elements ofp_(AE)(c)·⅓(p_(A,AE,min)(c)+p_(A,AE,max)(c)+p_(A,AE,mode)(c)),

/d/ a circuit for the determining of a value of expectation of height ofliquid hydrocarbon for said column analytically according to the sumdetermined by the circuit /c/ for each mesh c of said column and foreach facies A.

As such, this invention also aims for a computer program comprisinginstructions for the implementing of the method described hereinabove,when this program is executed by a processor.

This program can use any programming language (for example, an objectlanguage or other), and be in the form of interpretable source code, apartially compiled code or an entirely compiled code. FIG. 3 describedin detail hereinafter, can form the flow chart of the general algorithmof such a computer program.

BRIEF DESCRIPTION OF THE DRAWINGS

Other characteristics and advantages of the invention shall furtherappear when reading the following description. The latter is solely forthe purposes of illustration and must be read with respect to theannexed drawings wherein:

FIG. 1a shows a possible example of a geological model;

FIG. 1b describes a possible example of a probability distribution thathas a triangular distribution;

FIG. 2a describes an example of a transversal cross-section of areservoir comprising water, gas and oil as well as the rate of watersaturation of the rock of this reservoir.

FIG. 2b describes an example of a water saturation curve of the rock ofa reservoir modelled as being continuous and linear in pieces;

FIG. 3 is an example of a flow chart of an embodiment of the invention;

FIG. 4 shows an example of a device allowing for the implementation ofan embodiment of the invention;

FIGS. 5a and 5b are respectively examples of maps of expectation and ofvariance of height of liquid hydrocarbon in the model determined usingmulti-creation tools;

FIGS. 5c and 5d are respectively examples of maps of expectation and ofvariance of height of liquid hydrocarbon in the model determined usingmethods according to a possible embodiment of the invention.

DETAILED DESCRIPTION OF THE DRAWINGS

FIG. 1a shows a possible example of a transversal cross-section of ageological model 100.

This geological model 100 is comprised of a large number of pixels (ex.104, 105) each one describing a facies present in the geologicalsubsurface.

The petrol-physical data NTG, Phi and So are often determined for eachfacies. For each pixel c (of dimensions h_(p)xh_(p)xh_(p)), it ispossible to calculate a “pixel” value of HuPhiSo defined as follows:

$\sum\limits_{A_{i}}^{\;}\; {1_{\;_{A_{i}}}{\left( {x,y,z} \right) \cdot h_{p} \cdot \phi_{c,A_{i}} \cdot {NTG}_{c,A_{i}} \cdot {SO}_{c,A_{i}}}}$

where 1_(A) _(i) (x, y, z) is the signalling function of the faciesA_(i) (equal to 1 if the facies for the pixel of coordinates (x,y,z) isA_(i), equal to 0 otherwise).

Of course, each pixel is associated with a single facies and:

${\sum\limits_{A_{i}}^{\;}\; {1_{\;_{A_{i}}}\left( {x,y,z} \right)}} = 1$

Moreover, according to the geological history of the subsurface, thesedimentary filling and the geodynamic context that gave rise to therocks of the reservoir, the sedimentation can be done heterogeneouslybut in a structured manner. The sedimentary environments created as suchare more or less rich with certain facies.

For example, in a channel, the facies that are most present will be thecoarse sands. The farther the rock is from the channel the more clayeyit is.

These structures are generally called the architectural elements (AE).

For a given pixel of the model, it is considered that there is only onearchitectural element. The geological model 100 comprises threearchitectural elements 101, 102, and 103.

The interfaces between two architectural elements are often marked onmaps by geophysicists. For example, the interface 102 _(min) canrepresent a possible interface between the architectural element 102 andthe architectural element 103.

It is as such possible to introduce, as for the facies, a signalling ofarchitectural elements 1_(AE) _(j) (x, y, z) indicating the presence ornot of the architectural element AE_(j) for the pixel of coordinates(x,y,z).

The proportions of facies in a given pixel c are, most often, a functionof the architectural element associated with this pixel c. As such thesignalling of facies A_(i), for a given architectural element AE_(j), isexpressed in the form 1_(A) _(i) _(|AE) _(j) (x, y, z) (also noted as1_(A) _(i) _(,AE) _(j) (x, y, z) or 1_(A) _(i) _(,AE) _(j) (c)).

As such, the “pixel” HuPhiSo value is defined as follows:

$\sum\limits_{A_{i}}^{\;}{\sum\limits_{{AE}_{j}}^{\;}{\overset{\;}{1_{\;_{{AE}_{j}}}{\left( {x,y,z} \right) \cdot 1_{{Ai}{AE}_{j}}}}{\left( {x,y,z} \right) \cdot h_{c} \cdot \phi_{c,A_{i}} \cdot {NTG}_{c,A_{i}} \cdot {SO}_{c,A_{i}}}}}$

with h_(c) the height of the pixel c under consideration.

It is also possible to determine the “column” HuPhiSo values by addingthe “pixel” HuPhiSo values for all of the columns of the model (forexample, column 106). This “column” HuPhiSo value is given by thefollowing formula:

$\sum\limits_{c}^{\;}{\sum\limits_{A_{i}}^{\;}{\left( {\sum\limits_{{AE}_{j}}^{\;}{1_{{AE}_{j}}{(c) \cdot 1_{{Ai}{AE}_{j}}}(c)}} \right) \cdot h_{c} \cdot \phi_{c,A_{i}} \cdot {NTG}_{c,A_{i}} \cdot {SO}_{c,A_{i}}}}$

with c the pixels of the column under consideration and with h_(c) theheight of the current pixel.

In the general case, the signalling 1_(A) _(i) _(|AE) _(j) (x, y, z) isnot known in a precise manner, pixel by pixel: only the proportion cubesof the various facies A_(i) are known. In order to overcome thisproblem, prior solutions proposed:

-   -   to determine the limits of the architectural elements        stochastically;    -   to determine a map of facies by stochastically simulating a        possible value of facies for each pixel thanks to the given        proportion cubes for each facies of each architectural element,    -   to determine “pixel” HuPhiSo value of each pixel;    -   to determine “pixel HuPhiSo values for each column of the model        (i.e. HuPhiSo map);    -   to reiterate the first four steps a large number of times so as        to obtain an average and a converging variance of the HuPhiSo        maps.

This invention does not choose this stochastic path.

Indeed, it is possible to determine for each column of the model, theaverage of the column HuPhiSo values and the variance of the columnHuPhiSo values analytically, after certain hypotheses and geologicalsimplifications presented hereinafter.

Determination of the Expectation

The expectation of the column HuPhiSo values can be expressed in theform:

${{Esp}({HuPhiSO})} = {\sum\limits_{c = 1}^{n}{\sum\limits_{A_{i}}^{\;}{h_{c} \cdot \left( {\sum\limits_{{AE}_{j}}^{\;}\left( {{{Esp}\left\lbrack {1_{{AE}_{j}}(c)} \right\rbrack} \cdot {{Esp}\left\lbrack {1_{{Ai}{AE}_{j}}(c)} \right\rbrack}} \right)} \right) \cdot \left( {{{Esp}\left( \phi_{c,A_{i}} \right)} \cdot {{Esp}\left( {NTG}_{c,A_{i}} \right)} \cdot {{Esp}\left( {SO}_{c,A_{i}} \right)}} \right)}}}$

by assuming that the pixels c of the column are indexed from 1 to n(with n being a positive integer greater than 2).

The expression of this average creates certain geological hypotheses inorder to simplify the expression of this average: e.g. the independenceof the parameters φ_(c,A) _(i) , NTG_(c,A) _(i) and SO_(c,A) _(i) (assuch, the expectation of the product of these variables is also theproducts of their expectation). As such, this expression is not thesimple mathematical resolution of the problem posited.

If the signalling of facies 1_(A) _(i) (c) is a random variablefollowing a Bernouilli distribution of which the parameter is the localproportion of the facies, then the expectation of this variable knowingp_(A) _(i) (c) can be expressed in the following form:

Esp(1_(A) _(i) (c)|p _(A) _(i) (c))=Esp(1_(A) _(i) (c)=1|p _(A) _(i)(c))=p _(A) _(i) (c)

In conclusion:

Esp(1_(A) _(i) (c))=Esp(p _(A) _(i) (c))

The proportions p_(A) _(i) (c) are determined with certainty only onwells. There is therefore an uncertainty on the proportion of facies inthe rest of the reservoir. Because of this uncertainty, the proportionof facies can be, it too, a random variable, of which the probabilitydistribution can be modelled by a triangular distribution 150 describedin FIG. 1b (for a given facies A_(i) and for a pixel c) havingp_(min,Ai) as a minimum value, p_(mode,Ai) as a mode value andp_(max,Ai) as a maximum value. As such:

${{Esp}\left( {1_{{Ai}{AE}_{j}}(c)} \right)} = \frac{{p_{\min,{A_{i}{AE}_{j}}}(c)} + {p_{{mode},{A_{i}{AE}_{j}}}(c)} + {p_{\max,{A_{i}{AE}_{j}}}(c)}}{3}$

Of course, the reservoir is not necessarily stationary, i.e. itsproportions p_(A) _(i) (c) (and P_(min,A) _(i) (c),p_(mode,A) _(i)(c),p_(max,A) _(i) (c)) can vary according to the cell c underconsideration: for example, by approaching a limit or a particularpoint, a facies can tend to disappear or tend to predominate.

The same reasoning is valid for calculating the expectation of thesignalling of architectural elements. However, the proportion ofarchitectural elements does not have any uncertainty even if severalpositions of interfaces between architectural elements are possible (theoperator is often requested to point a minimum interface 102 _(min), anda maximum interface 102 _(max) (possibly a probable interface 102_(mode), if this latter interface is not determined, it is consideredthat this “probable interface” is located at an equal distance from theminimum and maximum interfaces). As such:

Esp(1_(AE) _(j) (c))=p _(AE) _(j) (c)

In conclusion, and through analogy, the expression Σ_(AE) _(j)Esp[1_(AE) _(j) (c)]. Esp1Ai|AEjc can be expressed in the form

$\Sigma_{{AE}_{j}}{{p_{{AE}_{j}}(c)} \cdot \frac{{p_{\min,{A_{i}{AE}_{j}}}(c)} + {p_{{mode},{A_{i}{AE}_{j}}}(c)} + {p_{\max,{A_{i}{AE}_{j}}}(c)}}{3}}$

which is in fact the proportion of facies A_(i) for the pixel c (whichtakes into account the possible presence of several AE and thereforeseveral proportion cubes for the same facies).

A rock is comprised of grains and of voids between the grains. Theporosity of the rock is the ratio between the volume of the voids andthe volume of the rock:

${\phi = \frac{V_{void}}{V_{total}}}\;$

This porosity is most often apportioned according to a naturaldistribution such as a normal distribution. Its average is, it too,uncertain and follows a triangular probability distribution.

Phi is therefore geologically the sum of two random variables:

φ=φ_(avg)+φ_(stoch)

where φ_(avg) is a random variable that can follow a triangulardistribution (similar to FIG. 1b ), of minimum value φ_(min), of modevalue φ_(mode) and of maximum value φ_(max), and where φ_(stoch) is arandom variable which follows a centred normal distribution (with zeroaverage).

As such:

${{Esp}(\phi)} = {{{{Esp}\left( \phi_{avg} \right)} + {{Esp}\left( \phi_{stoch} \right)}} = {{{Esp}\left( \phi_{avg} \right)} = \frac{\phi_{\min} + \phi_{mode} + \phi_{\max}}{3}}}$

The variable NTG represents the fraction of a rock volume favourable tothe presence of hydrocarbons. If V is a rock volume, then NTG*V is thevolume favourable to the presence of hydrocarbons, it is the usefulvolume of the rock. The value NTG can be apportioned according to anatural distribution that generally follows (but not necessarily) anormal distribution of average NTG_(avg).

The average NTG_(avg) can be, itself, uncertain and can follow atriangular probability distribution (similar to FIG. 1b ), of minimumvalue NTG_(min), of mode value NTG mode and of maximum value NTG_(max).

The distribution NTG is therefore geologically the sum of two randomvariables:

NTG=NTG_(avg)+NTG_(stoch)

where NTG_(avg) is a random variable which follows the triangulardistribution described hereinabove and where NTG_(stoch) is a randomvariable which follows a centred normal distribution (with zeroaverage).

${{Esp}({NTG})} = {{{Esp}\left( {NTG}_{avg} \right)} = \frac{{NTG}_{\min} + {NTG}_{mode} + {NTG}_{\max}}{3}}$

The variable SO represents the oil saturation of the rock.

The pores of the rocks can be filled with gas, water or oil (petrol).The oil saturation So, the water saturation Sw and the gas saturation Sgare therefore complementary:

So+Sw+Sg=1

In addition, the hydrocarbon saturation (i.e. gas and oil) is noted asS_(HC):

S _(HC) =So+S _(g)

There is therefore:

S _(HC)=1−Sw

In the reservoirs (porous and permeable rocks), the fluids are subjectedto gravity and apportioned according to their density. As such, under apermeable rock 201, it is possible to find (see FIG. 2a ) gas at the top(zone 202), oil/petrol in the middle (zone 203) and finally water (zone204).

Underneath the oil/water contact 205 (dimension z₂₀₅), the pores arefilled only with water and, above the contact 205 (zone 203), they canbe filled with oil but residual water can also be absorbed in the grainsof the rock.

As such, below the dimension z₂₀₅, the water saturation Sw is 1 (curve211). Above the dimension z₂₀₅, the water saturation Sw decreasesprogressively (curve 212) to a dimension (z_(R)) at which the watersaturation Sw remains substantially constant (i.e. the residual value ofwater saturation) and equal to V_(R) (curve 213).

The saturation curves (211, 212, 213) can be, themselves, uncertain (fora fixed dimension z₂₀₅ and known for certain). It is as such possible todetermine three water saturation curves:

-   -   a curve having a minimum water saturation (214),    -   a curve having a maximum water saturation (215),    -   a mode curve having the most likely case (212 and 213).

As such, for each height above the contact 205, it is possible to definea triangular probability distribution for the water saturation usingthree corresponding points of the curves defined as such.

Moreover, the dimension z_(w) of the interface between the oil and thewater (z₂₀₅ in the example of FIG. 2a ) may not be known for certain.This dimension z_(w) can also be represented by a random variable thatfollows a law of triangular probability of maximum value z_(w,min), ofminimum value z_(w,max) and of mode value Z_(w,mode).

S_(HC) (Z)=CS_(HC) (z−z_(w)) is noted for revealing the fact that thecurve of hydrocarbon saturation can be translated vertically due to theuncertainty on z_(w):

-   -   the curve having a minimum of hydrocarbon saturation is noted as        CS_(HC,min)(z−zw,    -   the curve having a maximum of hydrocarbon saturation is noted as        CS_(HC,max)(z−zw,    -   the mode curve having the most likely case is noted as        CS_(HC,mode)(z−z_(w)).

Likewise, the dimension z_(g) representing the dimension of the gas/oilinterface (z₂₀₆ in the example of FIG. 2a for the interface 206) cannotbe known for certain. This dimension z_(g) can be represented by arandom variable that follows a law of triangular probability of maximumvalue z_(g,min), of minimum value z_(g,max) and of mode value z_(g,mode)

As such, the expectation of the hydrocarbon saturation S_(HC) can beexpressed in the form:

$\begin{matrix}{{{Esp}\left( {S_{HC}(z)} \right)} = {{Esp}\left( {1 - {{Sw}(z)}} \right)}} \\{= {{Esp}\left( {{CS}_{HC}\left( {z - z_{w}} \right)} \right)}} \\{= {{Esp}\left( {{{Esp}\left( {{{CS}_{HC}\left( {z - z_{w}} \right)}z_{w}} \right)},z_{w}} \right)}}\end{matrix}$

As the variables CS_(HC) and Zw can be considered, geologically, as twoindependent random variables, it is possible to write:

Esp(S_(HC)(z)) = ∫_(z_(w, max ))^(z_(w, min ))∫_(cs_(HC, min ))^(cs_(HC, max ))cs_(HC)(z − z_(w)) ⋅ p(CS_(HC)(z − z_(w))z_(w)) ⋅ p(z_(w)) ⋅ dcs_(HC) ⋅ dz_(w)     Esp(S_(HC)(z)) = ∫_(z_(w, max ))^(z_(w, min ))p(z_(w)) ⋅ Esp(CS_(HC)(z − z_(w))  knowing  z_(w)) ⋅ dz_(w)${{Esp}\left( {S_{HC}(z)} \right)} = {\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{{cs}_{{HC},\min}\left( {z - z_{w}} \right)} + {{cs}_{{HC},{mode}}\left( {z - z_{w}} \right)} + {{cs}_{{HC},\max}\left( {z - z_{w}} \right)}}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}}$

As such, there are three terms to be calculated of the type

$\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}$

(with XXX among min, max and mode).

If the probability distribution of z_(w) is assumed to be triangular, itis possible to establish two cases:

if z _(w) ε[z _(w,max) ,z _(w,mode)]

${{p\left( z_{w} \right)} = {{{p_{1p} \cdot z_{w}} + p_{1c}} = {\frac{2 \cdot z_{w}}{\left( {z_{w,\min} - z_{w,\max}} \right) \cdot \left( {z_{w,{mode}} - z_{w,\max}} \right)}\; + \frac{{- 2}z_{w,\max}}{\left( {z_{w,\min} - z_{w,\max}} \right) \cdot \left( {z_{w,{mode}} - z_{w,\max}} \right)}}}}$if z _(w) ε[z _(w,mode) ,z _(w,min)]

${{p\left( z_{w} \right)} = {{{p_{2p} \cdot z_{w}} + p_{2c}} = {\frac{{- 2} \cdot z_{w}}{\left( {z_{w,\min} - z_{w,\max}} \right) \cdot \left( {z_{w,\min} - z_{w,{mode}}} \right)} + \mspace{11mu} \frac{2 \cdot z_{w,\min}}{\left( {z_{w,\min} - z_{w,\max}} \right) \cdot \left( {z_{w,\min} - z_{w,{mode}}} \right)}}}}\; $

As such:

${{\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}} = {{\int_{z_{w,\max}}^{z_{w,{mode}}}{\left\lbrack {\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}\left( {{p_{1p} \cdot z_{w}} + p_{1c}} \right)} \right\rbrack \cdot {dz}_{w}}}\; + {\int_{z_{w,{mode}}}^{z_{w,\min}}{\left\lbrack {\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}\left( {{p_{2p} \cdot z_{w}} + p_{2c}} \right)} \right\rbrack \cdot {dz}_{w}}}}}$

Then by carrying out a change in the variable h=z−z_(w), it is possibleto write:

${\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}} = {\int_{z - z_{w,{mode}}}^{z - z_{w,\max}}\left\lbrack {{\frac{{{cs}_{{HC},{XXX}}}^{(h)}}{3}{\left( {{p_{1p}\left( {z - h} \right)} + p_{1c}} \right) \cdot {dh}}} + {\int_{z - z_{w,\min}}^{z - z_{w,{mode}}}{\left\lbrack {\frac{{cs}_{{HC},{XXX}}(h)}{3}\left( {{p_{2p}\left( {z - h} \right)} + p_{2c}} \right)} \right\rbrack \cdot {dh}}}} \right.}$

In practice, the curve CS_(HC,XXX) (and therefore the curve Sw_(YYY)with YYY=min if XXX=max, YYY=mode if XXX=mode and YYY=max if XXX=min, asSw_(min=1)−CS_(HC,max)) is often linear by pieces, for example, betweentwo successive dimensions of the following lists:

-   -   h1₁, h1₂, . . . h1_(n-1), h1_(n1) (ordinates in increasing order        in the interval [z−z_(w,max),z−z_(w,mode)], with the end        dimensions corresponding to the boundaries of the interval, with        these end dimensions able to not correspond to a breakpoint of        the linearity of the curve, see hereinbelow) and    -   h2₁, h2₂, . . . h2_(n2-1), h2_(n2) (ordinates in increasing        order in the interval [z−z_(w,mode),z−zw,min, with the end        dimensions corresponding to the boundaries of the interval, as        such h2₁=h1_(n1), with these end dimensions able to not        correspond to a breakpoint of the linearity of the curve, see        hereinbelow).

For example, FIG. 2b shows an example of a water saturation curve1−CS_(HC,XXX)(z−z_(w)) (250), with the dimension z_(w) represented bythe line h=0.

This curve 250 is linear by piece between the successive points 251,252, 253, 254, 255, 256 and 257.

If the point of dimension (z does not correspond to one of the precedingpoints, it is possible to define a new point 271 on the linear portionof the water saturation curve between the points 252 and 253. The sameapplies for the dimension (z−z_(w,mode)) (corresponding then to the newpoint 272) and for the dimension (z−z_(w,min)) (corresponding then tothe new point 273).

As such:

${\int_{z - z_{w,\min}}^{z - z_{w,\max}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}} = {{\sum\limits_{i = 1}^{{n1} - 1}\; \left\lbrack {\int_{{h1}_{i + 1}}^{{h1}_{i}}{\left\lbrack {\frac{{cs}_{{HC},{XXX}}(h)}{3}\left( {{p_{1p}\left( {z - h} \right)} + p_{1c}} \right)} \right\rbrack \cdot {dh}}} \right\rbrack} + {\sum\limits_{i = 1}^{{n2} - 1}\left\lbrack {\int_{{h2}_{i + 1}}^{{h2}_{i}}{\left\lbrack {\frac{{cs}_{{HC},{XXX}}(h)}{3}\left( {{p_{2p}\left( {z - h} \right)} + p_{2c}} \right)} \right\rbrack \cdot {dh}}} \right\rbrack}}$

In the example of FIG. 2b , each one of the two preceding sums comprisestwo terms:

$\begin{matrix}{- {\sum\limits_{i = 1}^{{n1} - 1}\left\lbrack {\int_{{h1}_{i + 1}}^{{h1}_{i}}{\left\lbrack {\frac{{cs}_{{HC},{XXX}}(h)}{3}\left( {{p_{1p}\left( {z - h} \right)} + p_{1c}} \right)} \right\rbrack \cdot {dh}}} \right\rbrack}} & \;\end{matrix}$

(interval 281) comprises a terms corresponding to the integral betweenthe dimensions of the point 271 and of the point 253 and an integralterm between the dimensions of the point 253 and of the point 272;

$\begin{matrix}{- {\sum\limits_{i = 1}^{{n2} - 1}\left\lbrack {\int_{{h2}_{i + 1}}^{{h2}_{i}}{\left\lbrack {\frac{{cs}_{{HC},{XXX}}(h)}{3}\left( {{p_{2p}\left( {z - h} \right)} + p_{2c}} \right)} \right\rbrack \cdot {dh}}} \right\rbrack}} & \;\end{matrix}$

(interval 282) comprises a term corresponding to the integral betweenthe dimensions of the point 272 and of the point 254 and an integralterm between the dimensions of the point 254 and of the point 273.

It is considered that the equation describing the curve cs_(HC,XXX)between two points h1_(i) and h1_(i+1) is A1_(i,XXX)h+B1_(i,XXX) andthat the equation describing the curve cs_(HC,XXX) between two pointsh2_(i) and h2_(i+1) is A2_(i,XXX)h+B2_(i,XXXX).

In other words, the curve cs_(HC,XXX) is linear by piece between twofollowing successive points: (C1_(1,XXX), h1₁), (C1_(2,XXX), h1₂), . . .(C1_(n1-1,XXX), h1_(n1-1)), (C1_(n1,XXX), h1_(n1))=(C2_(1,XXX), h2₁),(C2_(2,XXX), h2₂), . . . (C2_(n2-1,XXX), h2_(n2-1)), (C2_(n2,XXX),h2_(n2)).

Then:

$\frac{{C1}_{{i + 1},{XXX}} - {C1}_{i,{XXX}}}{{h1}_{i + 1} - {h1}_{i}} = {{A1}_{{i,{XXX}}\mspace{14mu}}{and}}$$\mspace{14mu} {{{C1}_{i,{XXX}} - {\frac{{C1}_{{i + 1},{XXX}} - {C1}_{i,{XXX}}}{{h1}_{i + 1} - {h1}_{i}}{h1}_{i}}} = {{B1}_{i,{XXX}}{\mspace{11mu} \;}{and}}}$$\mspace{20mu} {\frac{{C2}_{{i + 1},{XXX}} - {C2}_{i,{XXX}}}{{h2}_{i + 1} - {h2}_{i}} = {{A2}_{{i,{XXX}}\mspace{14mu}}{and}}}\;$$\mspace{11mu} {{{C2}_{i,{XXX}} - {\frac{{C2}_{{i + 1},{XXX}} - {C2}_{i,{XXX}}}{{h2}_{i + 1} - {h2}_{i}}{h2}_{i}}} = {B2}_{i,{XXX}}}{\; \mspace{25mu}}$

and

As such:

$\begin{matrix}{{{\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}} = {{\sum\limits_{i = 1}^{{n1} - 1}\; \left\lbrack {\int_{{h1}_{i}}^{{h1}_{i + 1}}{\left\lbrack {\frac{{{A1}_{i,{XXX}} \cdot h} + {B1}_{i,{XXX}}}{3}\left( {{p_{1p}\left( {z - h} \right)} + p_{1c}} \right)} \right\rbrack \cdot {dh}}} \right\rbrack} + {\sum\limits_{i = 1}^{{n2} - 1}\left\lbrack {\int_{{h2}_{i}}^{{h2}_{i + 1}}{\left\lbrack {\frac{{{{A2}_{i}}_{,{XXX}} \cdot h} + {B2}_{i,{XXX}}}{3}\left( {{p_{2p}\left( {z - h} \right)} + p_{2c}} \right)} \right\rbrack \cdot {dh}}} \right\rbrack}}}\;} & \; \\\left. \left. {\left. {\left. \left. \mspace{79mu} {{{or}\text{:}}\mspace{79mu} {{\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}} = {{\frac{1}{3}{\sum\limits_{i = 1}^{{n1} - 1}\; \left\lbrack {{\left( {{{- p_{1p}} \cdot {A1}_{i,{XXX}} \cdot \frac{{h1}_{i + 1}^{3}}{3}} + \left( {{A1}_{i,{XXX}} \cdot \left( {p_{1c} + p_{1p}} \right) \cdot z} \right) - {p_{1p} \cdot {B1}_{i,{XXX}}}} \right) \cdot \frac{{h1}_{i + 1}^{2}}{2}} + {{B1}_{i,{XXX}} \cdot \left( {p_{1c} + {p_{1p}z}} \right) \cdot {h1}_{i + 1}}} \right)}} - {\left( {{{- p_{1p}} \cdot {A1}_{i,{XXX}} \cdot \frac{{h1}_{i}^{3}}{3}} + \left( {{{A1}_{i,{XXX}} \cdot \left( {p_{1c} + {p_{1p} \cdot z}} \right)} - {p_{1p} \cdot z}} \right) - {p_{1p} \cdot {B1}_{i,{XXX}}}} \right) \cdot \frac{{h1}_{i}^{2}}{2}} + {{B1}_{i,{XXX}} \cdot \left( {p_{1c} + {p_{1p}z}} \right) \cdot {h1}_{i}}}}} \right) \right\rbrack + {\frac{1}{3}{\sum\limits_{i = 1}^{{n2} - 1}\; {\left\lbrack {\left( {{{- p_{2p}} \cdot {A2}_{i,{XXX}} \cdot \frac{{h2}_{i + 1}^{3}}{3}} + {{A2}_{i,{XXX}} \cdot \left( {p_{2c} + p_{2p}} \right) \cdot z}} \right) - {p_{2p} \cdot {B2}_{i,{XXX}}}} \right) \cdot \frac{{h2}_{i + 1}^{2}}{2}}}} + {{B2}_{i,{XXX}} \cdot \left( {p_{2c} + {p_{2p}z}} \right) \cdot {h2}_{i + 1}}} \right) - {\left( {{{- p_{2p}} \cdot {A2}_{i,{XXX}} \cdot \frac{{h2}_{i}^{3}}{3}} + {{A2}_{i,{XXX}} \cdot \left( {p_{2c} + {p_{2p} \cdot z}} \right)} - {p_{2p} \cdot {B2}_{i,{XXX}}}} \right) \cdot \frac{{h2}_{i}^{2}}{2}} + {{B2}_{i,{XXX}} \cdot \left( {p_{2c} + {p_{2p}z}} \right) \cdot {h2}_{i}}} \right) \right\rbrack & \;\end{matrix}$

As such, each one of the terms of Esp(S_(HC)(z)) can be calculatedeasily

$\left( {i.e.{\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}}} \right.$

with XXX among min, max and mode).

Moreover, it is possible to define the “oil+water” signalling (or the“liquid-gas” signalling, with the liquid here being the oil and thewater) 1_(OW) by:

${1_{OW}(z)} = \left\{ \begin{matrix}{{1\mspace{14mu} {if}\mspace{14mu} z} < z_{g}} \\{0\mspace{14mu} {otherwise}}\end{matrix} \right.$

As such, the oil saturation can be expressed in the form:

So(z)=S _(HC)(z)·1_(OW)(z)

In conclusion:

Esp(So(z))=Esp(S _(HC)(z))·Esp(1_(OW)(z))

The “liquid-gas” signalling can be determined as follows:

${{Esp}\left( {1_{OW}(z)} \right)} = \left\{ \begin{matrix}{{0\mspace{14mu} {if}\mspace{14mu} z} \geq z_{g,\min}} \\{{1 - {\int_{Z_{G,\min}}^{z}{{p\left( z_{g} \right)}{dz}_{g}\mspace{14mu} {if}\mspace{14mu} z_{g,\max}}}} \leq z \leq z_{g,\min}} \\{{1\mspace{14mu} {if}\mspace{14mu} z} \leq z_{g,\max}}\end{matrix} \right.$

with p(z_(g)) the probability that the dimension of the gas-oilinterface z_(g) is located at this dimension, and with this probabilityhaving a distribution between the boundaries z_(g,max) and z_(g,min).

FIG. 3 shows a flow chart that allows for the determining of anexpectation of the HuPhiSo map, in an embodiment of the invention.

During the receiving of the receiving of the n_(x)m proportion cubes 301_(1,1), 301 _(2,1), . . . , 301 _(n,m) (respectively describing theproportions of the facies A₁|AE₁, A₂|AE₁, . . . A_(n)|AE_(m)) and of theapportionment cubes of the m architectural elements 302 ₁, 302 ₂, . . ., 302 _(m) (respectively describing the architectural element AE₁, AE₂,. . . AE_(m) in the model), it is possible to create (step 303) an emptyHuPhiSo model. This model is a model that comprises cells/pixels. Thismodel is of dimensions similar to the dimensions of the proportion cubesand/or of the apportionment cubes describing the apportionment of thearchitectural elements.

The proportion cube 301 _(i,j) comprises a plurality of cells c eachassociated with three proportion values as such describing a triangulardistribution of proportion: a minimum value p_(min,A) _(i) _(|AE) _(j)(c), a maximum value p_(max,A) _(i) _(|AE) _(j) (c) and a mode valuep_(mode,A) _(i) _(|AE) _(j) (c).

The apportionment cube 302 _(j) comprises a plurality of cells c witheach one associated with a proportion value p_(AE) _(j) (c).

If a cell of the model does not have any associated value of expectationof HuPhiSo (test 304, output OK), this cell is selected.

It is as such possible to determine (step 305) the expectation of theHuPhiSo value for the selected cell:

$\sum\limits_{A_{i}}^{\;}\; {h_{c} \cdot \left( {\sum\limits_{{AE}_{j}}^{\;}\; \left( {{p_{{AE}_{j}}(c)} \cdot \frac{{p_{\min,{A_{i}{AE}_{j}}}(c)} + {p_{{mode},{A_{i}{AE}_{j}}}(c)} + {p_{\max,{A_{i}{AE}_{j}}}(c)}}{3}} \right)} \right) \cdot \left( {{{Esp}\left( \phi_{c,A_{i}} \right)} \cdot {{Esp}\left( {NTG}_{c,A_{i}} \right)} \cdot {{Esp}\left( {SO}_{c,A_{i}} \right)}} \right)}$

This formula creates a certain number of simplifications based ongeological considerations mentioned hereinabove.

If it is considered that the variables φ_(c,A) _(i) , NTG_(c,A) _(i) ,or SO_(c,A) _(i) are certain (i.e. fixed values without uncertainty) forthe cells c and for the facies Ai, then the expectation of these valuesis equal to these certain values and the value determined in the step305 can be directly associated with the cell selected in the model.

If ever the variable φ_(c,A) _(i) is considered to be uncertain, threevalues that represent a triangular probability distribution are providedfor this cell selected and for each facies A_(i): a minimum valueφ_(c,A) _(i) _(,min), a maximum value φ_(c,A) _(i) _(,max), and a modevalue φ_(c,A) _(i) _(,mode). Then, the expectation of the porosity forthe selected cell and for the facies A_(i) is (step 306 a):

${{Esp}\left( \phi_{c,A_{i}} \right)} = \frac{\phi_{c,A_{i},\min} + \phi_{c,A_{i},{mode}} + \phi_{c,A_{i},\max}}{3}$

Likewise, if the variable NTG_(c,A) _(i) is considered to be uncertain,three values representing a triangular distribution of probability areprovided for this selected cell and for each facies A: a minimum valueNTG_(c,A) _(i) _(,min), a maximum value NTG_(c,A) _(i) _(,max), and amode value NTG_(c,A) _(i) _(,mode). Then, the expectation of the valueNTG for the selected cell and for the facies A_(i) is (step 306 b):

${{Esp}\left( {NTG}_{c,A_{i}} \right)} = \frac{{NTG}_{c,A_{i},\min} + {NTG}_{c,A_{i},{mode}} + {NTG}_{c,A_{i},\max}}{3}$

Finally, if the variable SO_(c,A) _(i) is considered to be uncertain,three water saturation curves are provided as a function of the depth zof the mesh c, a minimum, maximum and mode curve, linear by piece andrespectively defined by a plurality of points and by three valuesrepresenting a triangular distribution of probability for the value ofthe dimension of a water-oil interface: a minimum value z_(w,max), amaximum value z_(w,min), and a mode value Z_(w,mode). Then, theexpectation of the value SO for the selected cell is (step 306 c)calculated as indicated hereinabove.

A probability distribution p(z_(g)) of the oil-gas dimension z_(g) canalso be provided.

The value determined in the step 305 and in the step 306 a and/or 306 band/or 306 c can be used to calculate a value of the expectation of thevariable HuPhiSo. This expectation can then be associated with thiscell.

If all of the cells of the model have an associated value of expectationof HuPhiSo (test 304, output KO), it is possible to add (step 307) theexpectations of HuPhiSo of the cells that belong to the same column inorder to obtain an expectation of HuPhiSo column for said column underconsideration.

This map 308 is then provided to an operator for viewing and/or laterprocessing.

Determination of the Variance

Moreover, the variance of the HuPhiSo column values can be expressed inthe form:

Var(HuPhiSo_(column))=Esp(HuPhiSo_(column) ²)−Esp(HuPhiSo_(column))²

The calculation of the term Esp(HuPhiSo_(column)) and therefore ofEsp(HuPhiSo_(column))₂ is detailed hereinabove and the latter can becalculated for the calculation of the variance.

The formula hereinbelow makes the hypothesis that, geologically, thevariables φ_(c,A) _(i) , NTG_(c,A) _(i) and SO_(c,A) _(i) areindependent variables:

${{Esp}\left( {HuPhiSo}^{2} \right)} = {{\sum\limits_{A}^{\;}{\sum\limits_{B = A}^{\;}{\sum\limits_{c\; 1}^{\;}{\sum\limits_{{c\; 2} = {c\; 1}}^{\;}{h_{c\; 1}{h_{c\; 2} \cdot {{Esp}\left( {1_{A}\left( {c\; 1} \right)1_{B}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {\phi_{{c\; 1},A}\phi_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{NTG}_{{c\; 1},A}{NTG}_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{So}_{{c\; 1},A}{So}_{{c\; 2},B}} \right)}}}}}}} + {\sum\limits_{A}^{\;}{\sum\limits_{B \neq A}^{\;}{\sum\limits_{c\; 1}^{\;}{\sum\limits_{{c\; 2} = {c\; 1}}^{\;}{h_{c\; 1}{h_{c\; 2} \cdot {{Esp}\left( {1_{A}\left( {c\; 1} \right)1_{B}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {\phi_{{c\; 1},A}\phi_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{NTG}_{{c\; 1},A}{NTG}_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{So}_{{c\; 1},A}{So}_{{c\; 2},B}} \right)}}}}}}} + {\sum\limits_{A}^{\;}{\sum\limits_{B = A}^{\;}{\sum\limits_{c\; 1}^{\;}{\sum\limits_{{c\; 2} \neq {c\; 1}}^{\;}{h_{c\; 1}{h_{c\; 2} \cdot {{Esp}\left( {1_{A}\left( {c\; 1} \right)1_{B}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {\phi_{{c\; 1},A}\phi_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{NTG}_{{c\; 1},A}{NTG}_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{So}_{{c\; 1},A}{So}_{{c\; 2},B}} \right)}}}}}}} + {\sum\limits_{A}^{\;}{\sum\limits_{B \neq A}^{\;}{\sum\limits_{c\; 1}^{\;}{\sum\limits_{{c\; 2} \neq {c\; 1}}^{\;}{h_{c\; 1}{h_{c\; 2} \cdot {{Esp}\left( {1_{A}\left( {c\; 1} \right)1_{B}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {\phi_{{c\; 1},A}\phi_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{NTG}_{{c\; 1},A}{NTG}_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{So}_{{c\; 1},A}{So}_{{c\; 2},B}} \right)}}}}}}}}$

with A, B a facies (possibly identical) among the various faciespossible, and c1, c2 a cell (possible identical) of the column underconsideration for the calculation of the column HuPhiSO.

The first quadruple sum corresponds to the terms where the facies A andB are identical and the cells c1 and c2 are identical, the second sumcorresponds to the terms where the facies A and B are different but thecells c1 and c2 identical, the third sum corresponds to the cases wherethe facies A and B are the same but the cells c1 and c2 different andfinally the last sum corresponds to the cases where the facies A and Bare different and the cells c1 and c2 also.

The geological hypothesis is made that there cannot be two differentfacies in the same cell, i.e.:

Esp(1_(A)(c1)·1_(B)(c2))=0 if c1=c2 and A≠B

As such:

${{{Esp}\left( {HuPhiSo}^{2} \right)} = {{\sum\limits_{c}^{\;}{\sum\limits_{A}^{\;}{h_{c}^{2} \cdot {{Esp}\left( {1_{A}(c)^{2}} \right)} \cdot {{Esp}\left( \phi_{c,A}^{2} \right)} \cdot {{Esp}\left( {NTG}_{c,A}^{2} \right)} \cdot {{Esp}\left( {So}_{c,A}^{2} \right)}}}} + {\sum\limits_{A}^{\;}{\sum\limits_{c\; 1}^{\;}{\sum\limits_{{c\; 2} \neq {c\; 1}}^{\;}{h_{c\; 1}{h_{c\; 2} \cdot {{Esp}\left( {1_{A}\left( {c\; 1} \right)1_{A}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {\phi_{{c\; 1},A}\phi_{{c\; 2},A}} \right)} \cdot {{Esp}\left( {{NTG}_{{c\; 1},A}{NTG}_{{c\; 2},A}} \right)} \cdot {{Esp}\left( {{So}_{{c\; 1},A}{So}_{{c\; 2},A}} \right)}}}}}} + {\sum\limits_{A}^{\;}{\sum\limits_{B \neq A}^{\;}{\sum\limits_{c\; 1}^{\;}{\sum\limits_{{c\; 2} \neq {c\; 1}}^{\;}\; {h_{c\; 1}{h_{c\; 2} \cdot {{Esp}\left( {1_{A}\left( {c\; 1} \right)1_{B}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {\phi_{{c\; 1},A}\phi_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{NTG}_{c\; 1.A}{NTG}_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{So}_{{c\; 1},A}{So}_{{c\; 2},B}} \right)}}}}}}}}}\mspace{14mu}$

This sum reveals a certain number of identical terms that it is possibleto group together in order to decrease the required calculation time.

If the column under consideration comprises n cells and if there are Nfacies possible (the latter then being indexed), it is possible towrite:

$\left. {{{Esp}\left( {HuPhiSo}^{2} \right)} = {{\sum\limits_{c}^{\;}{\sum\limits_{A}^{\;}\left\lbrack {h_{c}^{2} \cdot {{Esp}\left( \phi_{c,A}^{2} \right)} \cdot {{Esp}\left( {NTG}_{c,A}^{2} \right)} \cdot {{Esp}\left( {So}_{c,A}^{2} \right)} \cdot {\sum\limits_{AE}^{\;}\; {{{Esp}\left( {p_{AE}(c)} \right)} \cdot {{Esp}\left( {p_{A,{AE}}(c)} \right)}}}} \right\rbrack}} + {2{\sum\limits_{A = 1}^{N}{\sum\limits_{{c\; 1} = 1}^{n - 1}{\sum\limits_{{c\; 2} = {{c\; 1} + 1}}^{n}\; \left\lbrack {h_{c\; 1}h_{c\; 2}{{{Esp}\left( {1_{A}\left( {c\; 1} \right)1_{A}\left( {c\; 2} \right)} \right)} \cdot {\quad{{Esp}\; {\left( {\phi_{{c\; 1},A}\phi_{{c\; 2},A}} \right) \cdot {{Esp}\left( {{NTG}_{{c\; 1},A}{NTG}_{{c\; 2},A}} \right)} \cdot}}\quad}}{{{Esp}\left( {{So}_{{c\; 1},A}{So}_{{c\; 2},A}} \right)} \cdot {\underset{{AE}\; 1}{\overset{\;}{\quad\sum}}{\sum\limits_{{AE}\; 2}^{\;}\; {{{Esp}\left( {1_{{AE}\; 1}\left( {c\; 1} \right)1_{{AE}\; 2}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {1_{A,{{AE}\; 1}}\left( {c\; 1} \right)1_{A,{{AE}\; 2}}({c2})} \right)}}}}}} \right\rbrack}}}} + {2{\sum\limits_{{c\; 1} = 1}^{n - 1}{\sum\limits_{{c\; 2} = {{c\; 1} + 1}}^{n}{\sum\limits_{A = 1}^{N}{\sum\limits_{B \neq A}^{N}{\left\lbrack {h_{c\; 1}h_{c\; 2}{{{Esp}\left( {1_{A}\left( {c\; 1} \right)1_{B}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {\phi_{{c\; 1},A}\phi_{{c\; 2},B}} \right)} \cdot}}\quad \right.{{{Esp}\left( {{NTG}_{{c\; 1},A}{NTG}_{{c\; 2},B}} \right)} \cdot {{Esp}\left( {{So}_{{c\; 1},A}{So}_{{c\; 2},B}} \right)} \cdot \; {\sum\limits_{{AE}\; 1}^{\;} {\sum\limits_{{AE}\; 2}^{\;}{{{Esp}\left( {1_{{AE}\; 1}\left( {c\; 1} \right) 1_{{AE}\; 2} \left( {c\; 2} \right)} \right)} \overset{\;}{\cdot \quad}{\quad\quad} {{Esp}\left( {1_{A,{{AE}\; 1}}\left( {c\; 1} \right)1_{B,{{AE}\; 2}}({c2})} \right)}}}}}}}}}}}}} \right\rbrack \mspace{40mu}$

The term Σ_(AE)Esp(p_(AE)(c))·Esp(p_(A,AE)(c)) was calculated previouslyand, following the geological simplifications made, is:

$\sum\limits_{AE}^{\;}\; \frac{{p_{AE}(c)} \cdot \left( {{p_{A,{{AE} \cdot \min}}(c)} + {p_{A,{{AE} \cdot {mode}}}(c)} + {p_{A,{{AE} \cdot \max}}(c)}} \right.}{3}$

The term Esp(1_(AE1)(c1)1_(AE2)(c2)) can be expressed as follows if thegeological hypothesis is made according to which the uncertaintiesconcerning the architectural elements are only horizontal (there are novertical uncertainties):

${{Esp}\left( {1_{{AE}\; 1}{\left( {c\; 1} \right) \cdot 1_{{AE}\; 2}}\left( {c\; 2} \right)} \right)} = \left\{ \begin{matrix}{{{p_{{AE}\; 1}\left( {c\; 1} \right)} \cdot {p_{{AE}\; 2}\left( {c\; 2} \right)}}\mspace{14mu} {if}\mspace{14mu} {at}\mspace{14mu} {least}\mspace{14mu} {one}{\mspace{11mu} \;}{of}\mspace{14mu} {the}\mspace{14mu} {AE}{\mspace{11mu} \;}{is}\mspace{14mu} {certain}} \\{{{p_{AE}\left( {c\; 1} \right)}\mspace{14mu} {if}\mspace{14mu} {the}{\mspace{11mu} \;}{AE}{\; \mspace{11mu}}{are}\mspace{14mu} {all}{\mspace{11mu} \;}{uncertain}{\mspace{11mu} \;}{and}\mspace{14mu} {AE}\; 1} = {{AE}\; 2}} \\{{0{\mspace{11mu} \;}{if}\mspace{14mu} {the}\mspace{14mu} {AE}\mspace{14mu} {are}\mspace{14mu} {all}\mspace{14mu} {uncertain}\mspace{14mu} {and}\mspace{14mu} {AE}\; 1} \neq {{AE}\; 2}}\end{matrix} \right.$

That is to say:

${{Esp}\left( {1_{{AE}\; 1}{\left( {c\; 1} \right) \cdot 1_{{AE}\; 2}}\left( {c\; 2} \right)} \right)} = \left\{ \begin{matrix}{{p_{{AE}\; 2}\left( {c\; 2} \right)}\mspace{14mu} {if}\mspace{14mu} \left( {{p_{{AE}\; 1}\left( {c\; 1} \right)} = 1} \right)} \\{{p_{{AE}\; 1}\left( {c\; 1} \right)}{\mspace{11mu} \;}{{if}{\mspace{11mu} \;}\left( {{p_{{AE}\; 2}\left( {c\; 1} \right)} = 1} \right)}} \\{0\mspace{14mu} {if}\mspace{14mu} \left( {{p_{{AE}\; 1}\left( {c\; 1} \right)} = 0} \right)\mspace{14mu} {or}{\mspace{14mu} \;}{if}\mspace{14mu} \left( {{p_{{AE}\; 2}\left( {c\; 2} \right)} = 0} \right)} \\{{otherwise}\mspace{14mu} \left\{ \begin{matrix}{{{p_{AE}\left( {c\; 1} \right)}\mspace{14mu} {if}\mspace{11mu} {AE}\; 1} = {{AE}\; 2}} \\{{0{\mspace{11mu} \;}{if}\mspace{11mu} {AE}\; 1} \neq {{AE}\; 2}}\end{matrix} \right.}\end{matrix} \right.$

Then, we have:

${\sum\limits_{{AE}\; 1}^{\;}{\sum\limits_{{AE}\; 2}^{\;}\; {{{Esp}\left( {1_{{AE}\; 1}{\left( {c\; 1} \right) \cdot 1_{{AE}\; 2}}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {1_{A,{{AE}\; 1}}{\left( {c\; 1} \right) \cdot 1_{B,{{AE}\; 2}}}\left( {c\; 2} \right)} \right)}}}} = {{\sum\limits_{{AE}\; 1}^{\;}\; {{{Esp}\left( {1_{{AE}\; 1}{\left( {c\; 1} \right) \cdot 1_{{AE}\; 1}}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {1_{A,{{AE}\; 1}}{\left( {c\; 1} \right) \cdot 1_{B,{{AE}\; 1}}}\left( {c\; 2} \right)} \right)}}} + {\sum\limits_{{AE}\; 1}^{\;}\; {\sum\limits_{{{AE}\; 2} \neq {{AE}\; 1}}^{\;}\; {{{Esp}\left( {1_{{AE}\; 1}{\left( {c\; 1} \right) \cdot 1_{{AE}\; 2}}\left( {c\; 2} \right)} \right)} \cdot {{Esp}\left( {1_{A,{{AE}\; 1}}{\left( {c\; 1} \right) \cdot 1_{B,{{AE}\; 2}}}\left( {c\; 2} \right)} \right)}}}}}$

It is assumed that, geologically, the signalling of facies of two cellslocated in different architectural elements are independent betweenthem. As such, for AE1≠AE2:

Esp(1_(A, AE 1)(c 1) ⋅ 1_(B, AE 2)(c 2)) = Esp(1_(A, AE 1)(c 1)) ⋅ Esp(1_(B, AE 2)(c 2)) = Esp(p_(A, AE 1)(c 1)) ⋅ Esp(p_(B, AE 2)(c 2))

Moreover, after certain geological considerations, for AE1=AE2, it ispossible to write, if the facies are identical (i.e. A=B):

${{Esp}\left( {1_{A,{AE}}{\left( {c\; 1} \right) \cdot 1_{A,{AE}}}\left( {c\; 2} \right)} \right)} = {{\left( {1 - {\rho \left( {A,{\Delta \; z}} \right)}} \right) \cdot \left( {{{{Esp}\left( p_{A,{AE},{c\; 1}} \right)} \cdot {{Esp}\left( p_{A,{AE},{c\; 2}} \right)}} + \sqrt{{{Var}\left( p_{A,{AE},{c\; 1}} \right)} \cdot {{Var}\left( p_{A,{AE},{c\; 2}} \right)}}} \right)} + {{\rho \left( {A,{\Delta \; z}} \right)} \cdot \left( {{\frac{1}{2}{{Esp}\left( p_{A,{AE},{c\; 1}} \right)}} + {\frac{1}{2}{{Esp}\left( p_{A,{AE},{c\; 2}} \right)}}} \right)}}$

with ρ(A,Δz) the “correlogram” function of the facies A (also noted asρ_(A)(Δz) or ρ_(A)(c1,c2) and defined hereinafter) and Δz the differenceof vertical dimensions according to the vertical axis z (i.e. the axisof the column) between the cells c₁ and c₂.

In addition, for different facies (i.e. A≠B), it is possible to write:

${{Esp}\left( {1_{A,{AE}}{\left( {c\; 1} \right) \cdot 1_{B,{AE}}}\left( {c\; 2} \right)} \right)} = {\left( {{{{Esp}\left( p_{A,{AE}} \right)} \cdot {{Esp}\left( p_{B,{AE}} \right)}} \mp \sqrt{{{Var}\left( {p_{A,{AE}}\left( {c\; 1} \right)} \right)} \cdot {{Var}\left( {p_{B,{AE}}\left( {c\; 2} \right)} \right)}}} \right) \cdot \left( {1 - {\rho \left( {A,{\Delta \; z}} \right)}} \right)}$

The term variogram is used to refer to a random variable X, the functionγ:

${\gamma \left( {X,h} \right)} = {\frac{1}{2} \cdot {{Esp}\left( \left( {{X\left( {u + h} \right)} - {X(u)}} \right)^{2} \right)}}$

This variogram γ represents the spatial structuring of the covariancebetween the variable X at a point of the space and the same variable Xat another point.

Most often, the variograms γ of a real-valued random variable respectthe following conditions:

-   -   The variogram γ is an increasing function of the distance, which        means that the more two points are separated then more the data        at these points is correlated.    -   Beyond a distance R_(γ) (also called “range”), the variogram γ        has reached its asymptotic value (or almost, i.e. 95% of this        value), which means that the spatial structuring of the data        exists only up to this characteristic distance R_(γ); beyond        this distance R_(γ), the data is no longer correlated.    -   For a zero distance, the variogram γ is zero, which means that        the correlation between the data is perfect for two points that        are very close to one another.

The correlogram ρ is called a standardised covariance:

${\rho \left( {X,h} \right)} = {\frac{{Cov}\left( {{X\left( {u + h} \right)},{X(u)}} \right)}{{Cov}\left( {{X(u)},{X(u)}} \right)} = {1 - \frac{\gamma \left( {X,h} \right)}{{Cov}\left( {{X(u)},{X(u)}} \right)}}}$

Most often, for a given real-valued random variable X (such as thefacies), the correlogram can be modelled in the form (with R_(γ) being apredetermined value):

-   -   of an exponential correlogram:

${\rho (h)} = {\exp \left( {- \frac{3h}{R_{\gamma}}} \right)}$

-   -   of a spherical correlogram:

${\rho (h)} = \left\{ \begin{matrix}{1 - {\frac{3}{2} \cdot \frac{h}{R_{\gamma}}} + {\frac{1}{2} \cdot \left( \frac{h}{R_{\gamma}} \right)^{3}}} & {{{if}\mspace{14mu} h} < R_{\gamma}} \\0 & {otherwise}\end{matrix} \right.$

-   -   of a Gaussian correlogram:

$\; {{\rho (h)} = \left\{ \begin{matrix}{\exp \left( {- \left( \frac{3h}{R_{\gamma}} \right)^{2}} \right)} & {{{if}\mspace{14mu} h} < R_{\gamma}} \\0 & {otherwise}\end{matrix} \right.}$

As mentioned hereinabove, geologically, the cells of the same facies canbe apportioned by blocs, which are often the result of the history ofthe sedimentary environment to which they belong.

For example, for a sedimentary environment of the channel type in themarine environment (sand, turbidities), the flows of sand can haveprovokes the formation of sand lobes. These lobes will be modelled bythe geologist as a group of cells of sand facies.

This trend that the cells of the same facies have to group themselvesinto groups, of the same shape and of the same average size, defines aspatial correlation of the facies. A cell that is close to another cellas such has more chance to be part of the same facies as the latter. Itis possible to model this spatial correlation with the variograms or acorrelogram:

Cov(1_(A)(c1),1_(A)(c2))=ρ_(A)(c1,c2)·Var(1_(A)(c1))·Var(1_(A)(c2))

The values R_(γ) or “ranges” (see hereinabove) of facies (in the threedimensions of space) are the dimensions of the geological forms observedin the field, in a given architectural element. A variogram or acorrelogram of the facies can then be defined by architectural element.

If the signalling of facies is correlated spatially via the variograms,the proportions of facies can be perfectly correlated between them.

For the calculation of the HuPhiSo, the variograms can be input data.Variograms are used in particular for spatially correlating thevariables NTG and Phi.

Geologically, the spatial correlations can apply only to the stochasticportion (noted as X_(stoch)) of the random variables (X), with theuncertainties on the average then being perfectly correlated betweenthem.

In addition, the correlations between the porosities or the variablesNTG for two different cells are carried out in the same facies. As such,there is a correlation between two porosities or variables NTG, in twodifferent cells, only for the same facies. For two different facies, thestochastic portions of the variables NTG and Phi can be independent.

It is then possible to calculate the covariances between the stochasticvariables for two different cells and the same facies:

Cov(NTG_(stoch)(c1),NTG_(stoch)(c2))

=ρ_(NTG)(c1,c2)·Var(NTG_(stoch)(c1))·Var(NTG_(stoch)(c2))

Cov(φ_(stoch)(c1),φ_(stoch)(c2))=ρ_(φ)(c1,c2)·Var(φ_(stoch)(c1))·Var(φ_(stoch)(c2))

or ρ_(NTG) is the vertical correlogram of the variable NTG and ρ_(φ) thevertical correlogram of the variable Phi.

With regards to the term Esp(φ_(c,A) ²), it is possible to write, aftercertain geological considerations, that Esp(φ_(c,A)²)=Esp(φ_(c,A))²+Var(φ_(c,A,avg))+Var(φ_(c,A,stoch)).

The same applies for the variable NTG and as such ESp(NTG_(c,A)²)=ESp(NTG_(c,A))²+Var(NTG_(c,A,moy))+Var(NTG_(c,A,stoch)).

With regards to the expectation of the product of the porosity for thesame facies and two different cells, it is possible to write, aftercertain geological considerations (in particular by assuming that thereis a perfect correlation between the uncertainties over the averageportion of the correlogram of the variable Phi):

Esp(φ_(c1,A)·φ_(c2,A))

=Esp(φ_(c1,A))·Esp(φ_(c2,A))+√{square root over(Var(φ_(c1,A,inc))·Var(φ_(c2,A,inc)))}

+ρ_(φ) _(stoch) √{square root over(Var(φ_(c1,A,stoch))·Var(φ_(c2,A,stoch)))}

or ρ_(φ) _(stoch) is the stochastic portion of the correlogram of thevariable Phi (see hereinabove) and φ_(c1,A,inc) is the uncertainty onthe average (i.e. if the distribution of φ_(c1,A) is triangular withminimum boundary a, maximum boundary c and mode b,φ_(c1,A,inc)=a²+b²+c²−ab−bc−ca/18).

The same applies for the variable NTG and as such:

Esp(NTG_(c1,A)·NTG_(c2,A))

=Esp(NTG_(c1,A))·Esp(NTG_(c2,A))+√{square root over(Var(NTG_(c1,A,inc))·Var(NTG_(c2,A,inc)))}

+ρ_(NTG) _(stoch) √{square root over(Var(NTG_(c1,A,stoch))·Var(NTG_(c2,A,stoch)))}

where ρ_(NTG) _(stoch) is the stochastic portion of the correlogram ofthe variable NTG (see hereinabove).

Likewise, if A≠B, it is then possible to write:

${{Esp}\left( {\phi_{{c\; 1},A} \cdot \phi_{{c\; 2},B}} \right)} = {{{{Esp}\left( \phi_{{c\; 1},A} \right)} \cdot {{Esp}\left( \phi_{{c\; 2},B} \right)}} + \sqrt{{{Var}\left( \phi_{{c\; 1},A,{inc}} \right)} \cdot {{Var}\left( \phi_{{c\; 2},B,{inc}} \right)}}}$${{Esp}\left( {{NTG}_{{c\; 1},A} \cdot {NTG}_{{c\; 2},B}} \right)} = {{{{Esp}\left( {NTG}_{{c\; 1},A} \right)} \cdot {{Esp}\left( {NTG}_{{c\; 2},B} \right)}} + \sqrt{{{Var}\left( {NTG}_{{c\; 1},A,{inc}} \right)} \cdot {{Var}\left( {NTG}_{{c\; 2},B,{inc}} \right)}}}$

As such, by re-using the notations hereinabove:

${{Esp}\left( {S_{HC}(z)}^{2} \right)} = {\frac{1}{6}\left( {{{Esp}\left( {S_{{HC},\min}(z)}^{2} \right)} + {{Esp}\left( {S_{{HC},{mode}}(z)}^{2} \right)} + {{Esp}\left( {S_{{HC},\max}(z)}^{2} \right)} + \left( {{{Esp}\left( {{S_{{HC},\min}(z)} \cdot {S_{{HC},{mode}}(z)}} \right)} + {{Esp}\left( {{S_{{HC},\min}(z)} \cdot {S_{{HC},\max}(z)}} \right)} + {{Esp}\left( {{S_{{HC},{mode}}(z)} \cdot {S_{{HC},\max}(z)}} \right)}} \right)} \right.}$

As an illustration, the value Esp (S_(HC,XX1)(z)·S_(HC,XX2)(z)) can becalculated (with XX1 and XX2 among mode, min and max and with XX1 andXX2 different) as follows:

${{Esp}\left( {{S_{{HC},{{XX}\; 1}}(z)} \cdot {S_{{HC},{{XX}\; 2}}(z)}} \right)} = {\sum\limits_{i = 1}^{{n\; 1} - 1}\left\{ {\left\lbrack {{{{- p_{1p}} \cdot A}\; {1_{i,{{XX}\; 1}} \cdot A}\; {1_{i,{{XX}\; 2}} \cdot \frac{h\; 1_{i + 1}^{4}}{4}}} + {\left\{ {{A\; {1_{i,{{XX}\; 1}} \cdot A}\; {1_{i,{{XX}\; 2}} \cdot \left( {{p_{1p} \cdot z} + p_{1c}} \right)}} - {p_{1p} \cdot \left( {{A\; 1_{i,{{XX}\; 1}}B\; 1_{i,{{XX}\; 2}}} + {A\; 1_{i,{{XX}\; 2}}B\; 1_{i,{{XX}\; 1}}}} \right)}} \right\} \cdot \frac{h\; 1_{i + 1}^{3}}{3}} + {\left\{ {\left( {{A\; 1_{i,{{XX}\; 1}}B\; 1_{i,{{XX}\; 2}}} + {A\; 1_{i,{{XX}\; 2}}B\; 1_{i,{{XX}\; 1}}}} \right) \cdot \left( {{p_{1p} \cdot z} + p_{1c}} \right)} \right\} \cdot \frac{h\; 1_{i + 1}^{2}}{2}} - {{p_{1c} \cdot B}\; {1_{i,{{XX}\; 1}} \cdot B}\; {1_{i,{{XX}\; 2}} \cdot h}\; 1_{i + 1}}} \right\rbrack - \left. \quad\left\lbrack {{{{- p_{1p}} \cdot A}\; {1_{i,{{XX}\; 1}} \cdot A}\; {1_{i,{{XX}\; 2}} \cdot \frac{h\; 1_{i}^{4}}{4}}} + {\left\{ {{A\; {1_{i,{{XX}\; 1}} \cdot A}\; {1_{i,{{XX}\; 2}} \cdot \left( {{p_{1p} \cdot z} + p_{1c}} \right)}} - {p_{1p} \cdot \left( {{A\; 1_{i,{{XX}\; 1}}B\; 1_{i,{{XX}\; 2}}} + {A\; 1_{i,{{XX}\; 2}}B\; 1_{i,{{XX}\; 1}}}} \right)}} \right\} \cdot \frac{h\; 1_{i}^{3}}{3}} + {\left\{ {\left( {{A\; 1_{i,{{XX}\; 1}}B_{1_{i,{{XX}\; 2}}}} + {A\; 1_{i,{{XX}\; 2}}B\; 1_{i,{{XX}\; 1}}}} \right) \cdot \left( {{p_{1p} \cdot z} + p_{1c}} \right)} \right\} \cdot \frac{h\; 1_{i}^{2}}{2}} - {{p_{1c} \cdot B}\; {1_{i,{{XX}\; 1}} \cdot B}\; {1_{i,{{XX}\; 2}} \cdot h}\; 1_{i}}} \right\rbrack \right\} + {\sum\limits_{i = 1}^{{n\; 2} - 1}\left\{ {\left\lbrack {{{{- p_{2p}} \cdot A}\; {2_{i,{{XX}\; 1}} \cdot A}\; {2_{i,{{XX}\; 2}} \cdot \frac{h\; 2_{i + 1}^{4}}{4}}} + {\left\{ {{A\; {2_{i,{{XX}\; 1}} \cdot A}\; {2_{i,{{XX}\; 2}} \cdot \left( {{p_{2p} \cdot z} + p_{2c}} \right)}} - {p_{2p} \cdot \left( {{A\; 2_{i,{{XX}\; 1}}B\; 2_{i,{{XX}\; 2}}} + {A\; 2_{i,{{XX}\; 2}}B\; 2_{i,{{XX}\; 1}}}} \right)}} \right\} \cdot \frac{h\; 2_{i + 1}^{3}}{3}} + {\left\{ {\left( {{A\; 2_{i,{{XX}\; 1}}B\; 2_{i,{{XX}\; 2}}} + {A\; 2_{i,{{XX}\; 2}}B\; 2_{i,{{XX}\; 1}}}} \right) \cdot \left( {{p_{2p} \cdot z} + p_{2\; c}} \right)} \right\} \cdot \frac{h\; 2_{i + 1}^{2}}{2}} - {{p_{2c} \cdot B}\; {2_{i,{{XX}\; 1}} \cdot B}\; {2_{i,{{XX}\; 2}} \cdot h}\; 2_{i + 1}}} \right\rbrack - {\quad{\quad\left\lbrack {{{- p_{2p}} \cdot A}\; {2_{i,{{XX}\; 1}} \cdot \left. \quad{{A {2_{i,{{XX}\; 2}} \cdot \frac{h\; 2_{i}^{4}}{4}}} + {\left\{ {{A\; {2_{i,{{XX}\; 1}} \cdot A}\; {2_{i,{{XX}\; 2}} \cdot \left( {{p_{2p} \cdot z} + p_{2c}} \right)}} - {p_{2p} \cdot \left( {{A\; 2_{i,{{XX}\; 1}}B\; 2_{i,{{XX}\; 2}}} + {A\; 2_{i,{{XX}\; 2}}B\; 2_{i,{{XX}\; 1}}}} \right)}} \right\} \cdot \frac{h\; 2_{i}^{3}}{3}} + {\left\{ {\left( {{A\; 2_{i,{{XX}\; 1}}B\; 2_{i,{{XX}\; 2}}} + {A\; 2_{i,{{XX}\; 2}}B\; 2_{i,{{XX}\; 1}}}} \right) \cdot \left( {{p_{2p} \cdot z} + p_{2c}} \right)} \right\} \cdot \frac{h\; 2_{i}^{2}}{2}} - {{p_{2c} \cdot B}\; {2_{i,{{XX}\; 1}} \cdot B}\; {2_{i,{{XX}\; 2}} \cdot h}\; 2_{i}}} \right\rbrack}} \right\}}}} \right.}} \right.}$

Moreover, with z₁ different from z₂:

${{Esp}\left( {{S_{{HC},A}\left( z_{1} \right)} \cdot {S_{{HC},B}\left( z_{2} \right)}} \right)} = {{\frac{1}{6}\left( {{E\left( {{S_{{HC},A,\min}\left( z_{1} \right)} \cdot {S_{{HC},B,\min}\left( z_{2} \right)}} \right)} + {E\left( {{S_{{HC},A,{mode}}\left( z_{1} \right)} \cdot {S_{{HC},B,{mode}}\left( z_{2} \right)}} \right)} + {E\left( {{S_{{HC},A,\max}\left( z_{1} \right)} \cdot {S_{{HC},B,\max}\left( z_{2} \right)}} \right)}} \right)} + {\frac{1}{12}\left( {{E\left( {{S_{{HC},A,\min}\left( z_{1} \right)} \cdot {S_{{HC},B,{mode}}\left( z_{2} \right)}} \right)} + {E\left( {{S_{{HC},A,{mode}}\left( z_{1} \right)} \cdot {S_{{HC},B,\min}\left( z_{2} \right)}} \right)} + {E\left( {{S_{{HC},A,\min}\left( z_{1} \right)} \cdot {S_{{HC},B,\max}\left( z_{2} \right)}} \right)} + {E\left( {{S_{{HC},A,\max}\left( z_{1} \right)} \cdot {S_{{HC},B,\min}\left( z_{2} \right)}} \right)} + {E\left( {{S_{{HC},A,{mode}}\left( z_{1} \right)} \cdot {S_{{HC},B,\max}\left( z_{2} \right)}} \right)} + {E\left( {{S_{{HC},A,\max}\left( z_{1} \right)} \cdot {S_{{HC},B,{mode}}\left( z_{2} \right)}} \right)}} \right)}}$

Each one of the terms of this formula can be calculated as follows (withXX1 and XX2 among mode, min and max):

${E\left( {{S_{{HC},A,{{XX}\; 1}}\left( z_{1} \right)} \cdot {S_{{HC},B,{{XX}\; 2}}\left( z_{2} \right)}} \right)} = {{\sum\limits_{i = 1}^{n\; 1}{\int_{z_{w,i}}^{z_{w,{i + 1}}}{\left( {{p_{1p} \cdot z_{w}} + p_{1c}} \right) \cdot \left( {{A\; {1_{i,{{XX}\; 1}} \cdot \left( {z_{1} - z_{w}} \right)}} + {B\; 1_{i,{{XX}\; 1}}}} \right) \cdot \left( {{A\; {2_{i,{{XX}\; 2}} \cdot \left( {z_{2} - z_{w}} \right)}} + {B\; 2_{i,{{XX}\; 2}}}} \right) \cdot {dz}_{w}}}} + {\sum\limits_{i = 2}^{n\; 2}{\int_{z_{w,i}}^{z_{w,{i + 1}}}{\left( {{p_{2\; p} \cdot z_{w}} + p_{2c}} \right) \cdot \left( {{A\; {1_{i,{{XX}\; 1}} \cdot \left( {z_{1} - z_{w}} \right)}} + {B\; 1_{i,{{XX}\; 1}}}} \right) \cdot \left( {{A\; {2_{i,{{XX}\; 2}} \cdot \left( {z_{2} - z_{w}} \right)}} + {B\; 2_{i,{{XX}\; 2}}}} \right) \cdot {dz}_{w}}}}}$

In addition:

Esp(So _(c,A))=Esp(S _(HC,A)(c))·Esp(1_(OW)(c))

Esp(So _(c,A) ²)=Esp(S _(HC,A)(c)²)·Esp(1_(OW)(c))

and

Esp(So _(c1,A) ·So _(c2,B))=Esp(S _(HC,A)(c1)·S_(HC,B)(c2))·Esp(1_(OW)(c1)·1_(OW)(c2))

knowing that:

Esp(1_(OW)(z ₁)·1_(OW)(z ₂))=Esp(1_(OW)(Max(z ₁ ,z ₂)))

Furthermore, there may be certain uncertainties as to the position ofthe structures in the reservoir, which as such can modify the “HuPhiSo”maps (of variance or of expectation).

These uncertainties can be linked

-   -   to the errors and to the inaccuracy of the measurements, for        example concerning the position and the orientation of surfaces.        This uncertainty primarily concerns the geometry of the        structures;    -   to the natural variability of the geological objects, which is        more or less substantial according to the geological objects        (e.g. the natural variability of a fault surface is much lower        than that of a mineralisation surface of a uranium deposit);    -   to the lack of knowledge, on the existence of a structure for        example.

It is possible to account for these uncertainties by integrating theminto the uncertainty on the position of the oil-water interfacementioned hereinabove.

By using the preceding notations, it is possible to consider that theuncertainty of the oil-water interface follows a triangular distribution(uncertainty mentioned hereinabove) but to which is added an uncertaintyof the normal distribution of deviation σ(z) (linked to the uncertaintyof the structure).

Under geological hypotheses, the term:

$\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}$

(mentioned hereinabove) can be replaced with, a change in the variableh=z−z_(w) having been carried out:

${\int_{z - {3\; {\sigma {(z)}}}}^{z + {3\; {\sigma {(z)}}}}{\int_{z - z_{w,{mode}}}^{z - z_{w,\max}}{{\left( {{p_{1p} \cdot \left( {z^{\prime} - h} \right)} + p_{1c}} \right) \cdot \frac{{cs}_{{HC},{XXX}}(h)}{3} \cdot {dh} \cdot \left( \frac{e^{{- \frac{1}{2}} \cdot {(\frac{z^{\prime}}{\sigma {(z)}})}^{2}}}{2\; \pi \; {\sigma (z)}} \right)}{dz}^{\prime}}}} + {\int_{z - {3\; {\sigma {(z)}}}}^{z + {3\; {\sigma {(z)}}}}{\int_{z - z_{w,\min}}^{z - z_{w,{mode}}}{{\left( {{p_{2p} \cdot \left( {z^{\prime} - h} \right)} + p_{2c}} \right) \cdot \frac{{cs}_{{HC},{XXX}}(h)}{3} \cdot {dh} \cdot \left( \frac{e^{{- \frac{1}{2}} \cdot {(\frac{z^{\prime}}{\sigma {(z)}})}^{2}}}{2\; \pi \; {\sigma (z)}} \right)}{dz}^{\prime}}}}$

As such, by using the preceding notations, it is possible to replace:

$\mspace{20mu} {\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{cs}_{{HC},{XXX}}\left( {z - z_{w}} \right)}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}}$  with:$\int_{z - {3\; {\sigma {(z)}}}}^{z + {3\; {\sigma {(z)}}}}{\left\{ {{\frac{1}{3}{\sum\limits_{i = 1}^{{n\; 1} - 1}\left\lbrack {\left( {{{{- p_{1p}} \cdot A}\; {1_{i,{XXX}} \cdot \frac{h\; 1_{i + 1}^{3}}{3}}} + {\left( {{A\; {1_{i,{XXX}} \cdot \left( {p_{1c} + {p_{1p} \cdot z^{\prime}}} \right)}} - {{p_{1p} \cdot B}\; 1_{i,{XXX}}}} \right) \cdot \frac{h\; 1_{i + 1}^{2}}{2}} + {B\; {1_{i,{XXX}} \cdot \left( {p_{1c} + {p_{1p}z^{\prime}}} \right) \cdot h}\; 1_{i + 1}}} \right) - \left( {{{{- p_{1p}} \cdot A}\; {1_{i,{XXX}} \cdot \frac{h\; 1_{i}^{3}}{3}}} + {\left( {{A\; {1_{i,{XXX}} \cdot \left( {p_{1c} + {p_{1p} \cdot z^{\prime}}} \right)}} - {{p_{1p} \cdot B}\; 1_{i,{XXX}}}} \right) \cdot \frac{h\; 1_{i}^{2}}{2}} + {B\; {1_{i,{XXX}} \cdot \left( {p_{1c} + {p_{1p}z^{\prime}}} \right) \cdot h}\; 1_{i}}} \right)} \right\rbrack}} + {\frac{1}{3}{\sum\limits_{i = 1}^{{n\; 2} - 1}\left\lbrack {\left( {{{{- p_{2p}} \cdot A}\; {2_{i,{XXX}} \cdot \frac{h\; 2_{i + 1}^{3}}{3}}} + {\left( {{A\; {2_{i,{XXX}} \cdot \left( {p_{2c} + {p_{2p} \cdot z^{\prime}}} \right)}} - {{p_{2p} \cdot B}\; 2_{i,{XXX}}}} \right) \cdot \frac{h\; 2_{i + 1}^{2}}{2}} + {B\; {2_{i,{XXX}} \cdot \left( {p_{2c} + {p_{2p}z^{\prime}}} \right) \cdot h}\; 2_{i + 1}}} \right) - \left( {{{{- p_{2p}} \cdot A}\; {2_{i,{XXX}} \cdot \frac{h\; 2_{i}^{3}}{3}}} + {\left( {{A\; {2_{i,{XXX}} \cdot \left( {p_{2c} + {p_{2p} \cdot z^{\prime}}} \right)}} - {{p_{2p} \cdot B}\; 2_{i,{XXX}}}} \right) \cdot \frac{h\; 2_{i}^{2}}{2}} + {B\; {2_{i,{XXX}} \cdot \left( {p_{2c} + {p_{2p}z^{\prime}}} \right) \cdot h}\; 2_{i}}} \right)} \right\rbrack}}} \right\} \left( \frac{e^{{- \frac{1}{2}} \cdot {(\frac{z^{\prime}}{\sigma {(z)}})}^{2}}}{2\; \pi \; {\sigma (z)}} \right){dz}^{\prime}}$

Moreover, it is possible to calculate the variance of the value HuPhiSotaking account of the uncertainties on the structure (orVar(HuPhiSo_(str) _(_) _(inc)), using the value of the variance of thisvalue HuPhiSo without uncertainty (or Var(HuPhiSo) calculatedhereinabove):

${{Var}\left( {HuPhiSo}_{str\_ inc} \right)} = {{{Var}({HuPhiSo})} + {\left\lbrack \frac{{Esp}\left( {HuPhiSo}_{str\_ inc} \right)}{{Esp}\left( H_{str\_ inc} \right)} \right\rbrack^{2}{{Var}\left( H_{str\_ inc} \right)}}}$

with H_(str) _(_) _(inc) the height of the column above the water-oilinterface.

FIG. 4 shows an example of a device for determining a map of expectationof height of liquid hydrocarbon in an embodiment of the invention.

In this embodiment, the device comprises a computer 400, comprising amemory 405 for storing instructions allowing for the implementation ofthe method, the data of the measurements received, and temporary datafor carrying out the various steps of the method such as describedhereinabove.

The computer further comprises a circuit 404. This circuit can be, forexample:

-   -   a processor able to interpret instructions in the form of a        computer program, or    -   an electronic card of which the steps of the method of the        invention are described in the silicon, or    -   a programmable electronic chip such as a FPGA chip        (“Field-Programmable Gate Array”).

This computer comprises an input interface 403 for the receiving of datain particular the proportion cubes required for the calculations of theinvention, and an output interface 406 for the supplying of maps ofexpectation of the height of liquid hydrocarbon. Finally, the computercan comprise, in order to allow for an easy interaction with a user, ascreen 401 and a keyboard 402. Of course, the keyboard is optional, inparticular in the framework of a computer having the form of a touchtablet, for example.

Moreover, the block diagram shown in FIG. 3 is a typical example of aprogram of which certain instructions can be carried out with theequipment described. As such, FIG. 3 can correspond to the flow chart ofthe general algorithm of a computer program in terms of the invention.

FIG. 5a is an example of a map of the expectation of the variableHuPhiSo for a meshed model calculated using methods of prior art known(i.e. calculation using a plurality of creations or “multi-creations”,with the expectation of these creations then being calculated).

FIG. 5b is an example of a map of the variance of the variable HuPhiSofor a meshed model calculated using methods of prior art known(multi-creations).

The meshed model, considered in these examples, comprises 80,000 meshes.The uncertainty concerning in particular the proportion cubes, theporosity, the water saturations are taken into account. The calculationof a hundred creations takes about 20h. However, it is usual todetermine rather 2000 creations during such a calculation as“multi-creation” in order to ensure a correct convergence of theresults. As such, the calculation time can largely exceed 20h.

For this same model and the same taking into account of uncertainties,the FIG. 5c is an example of a map of the expectation of the variableHuPhiSo determined using a method that implements one of the embodimentsof the invention. FIG. 5d is an example of a map of the variance of thevariable HuPhiSo determined using a method that implements one of theembodiments of the invention.

The calculation of the maps 5 c and 5 d using one of the embodiments ofthe invention takes about 2h for a result that is entirely comparablewith the maps 5 a and 5 b (respectively 2% difference on the average and0.5% differences on the average).

Of course, this invention is not limited to the forms of creationdescribed hereinabove by way of examples; it extends to otheralternatives.

Other embodiments are possible.

1. A method for determining a map of expectation of height of liquidhydrocarbon in a geological model, the geological model comprisingmeshes that can be associated with a facies among a plurality of faciesand with an architectural element among a plurality of architecturalelements, the method comprises, for at least one column of thegeological model: /a/ receiving an apportionment cube for eacharchitectural element AE among the plurality of architectural elements,with the apportionment cube describing for each mesh c of the model theprobability p_(AE)(c) that said mesh belongs to said architecturalelement, /b/ receiving a proportion cube for each facies among theplurality of facies and for each architectural element among theplurality of architectural elements, with the proportion cube associatedwith a facies A and with an architectural element AE describing for eachmesh c of the model the proportion of said facies in said mesh if saidmesh belongs to said architectural element AE, with said proportionbeing a random variable of triangular distribution defined by a minimumvalue p_(A,AE,min)(c), a maximum value p_(A,AE,max)(c), and a modep_(A,AE,mode)(c), /c/ determining, for each mesh c of said column andfor each facies A of the plurality of facies, the sum over the pluralityof architectural elements ofp_(AE)(c)·⅓(p_(A,AE,min)(c)+pA,AE,maxc+pA,AE,modec, and /d/ determininga value of expectation of height of liquid hydrocarbon for said columnanalytically according to the sum determined in the step /c/ for eachmesh c of said column and for each facies A.
 2. The method according toclaim 1, wherein, with each mesh c of the model being associated with aproportion of porosity for each one of the facies of said plurality offacies, with said proportion of porosity being a random variable oftriangular distribution defined by a minimum value φ_(A,min)(c), amaximum value φ_(A,max)(c), and a mode φ_(A,mode)(c), the method furthercomprises: /c-1/ determining, for each mesh c of said column and foreach facies A among the plurality of facies, the value⅓(φ_(A,min)(c)+φ_(A,max)(c)+φ_(A,mode)(c)), wherein, the determinationof the step /d/ is a function of the value determined in the step/c-1/for each mesh c of said column and for each facies A among theplurality of facies.
 3. The method according to claim 1, wherein, witheach mesh c of the model being associated with a proportion of rockvolume favourable to the presence of hydrocarbons for each one of thefacies of said plurality of facies, with said proportion of rock volumefavourable to the presence of hydrocarbons being a random variable oftriangular distribution defined by a minimum value NTG_(A,min)(c), amaximum value NTG_(A,max)(c), and a mode NTG_(A,mode)(c), the methodfurther comprises: /c-2/ determining, for each mesh c of said column andfor each facies A among the plurality of facies, the value⅓(NTG_(A,min)(c)+NTG_(A,max)(c)+NTG_(A,mode)(c)), wherein, thedetermination of the step /d/ is a function of the value determined inthe step /c-2/ for each mesh c of said column and for each facies Aamong the plurality of facies.
 4. The method according to claim 1,wherein, with each mesh c of the model being associated with a liquidhydrocarbon saturation, said liquid hydrocarbon saturation beingcharacterised by a liquid hydrocarbon-water interface dimension z_(w)having a distribution probability p(z_(w)) between a maximum valuez_(w,min), and a minimum value z_(w,max), and by a triangulardistribution of the hydrocarbon saturation for each dimension h abovethis interface dimension, said distribution of hydrocarbon saturationbeing defined by a minimum value CS_(HC,A,min)(h), a maximum valueCS_(HC,A,max)(h), and a mode CS_(HC,A,mode)(h), the method furthercomprises: /c-3/ determining, for each mesh c of said column having adimension z and for each facies A among the plurality of facies, thevalue$\int_{z_{w,\max}}^{z_{w,\min}}{\frac{{{cs}_{{HC},\min}\left( {z - z_{w}} \right)} + {{cs}_{{HC},{mode}}\left( {z - z_{w}} \right)} + {{cs}_{{HC},\max}\left( {z - z_{w}} \right)}}{3}{{p\left( z_{w} \right)} \cdot {dz}_{w}}}$wherein, the determination of the step /d/ is a function of the valuedetermined in the step /c-3/ for each mesh c of said column and for eachfacies A among the plurality of facies.
 5. The method according to claim4, wherein, with said liquid hydrocarbon saturation being characterisedby a liquid hydrocarbon-gas interface dimension z_(g) having aprobability p(z_(g)), of distribution between a maximum value z_(g,min),and a minimum value z_(g,max), the method further comprises: /c-4/determining, for each mesh c of said column having a dimension z, thevalue $\left\{ {\begin{matrix}0 & {{{if}\mspace{14mu} z} \geq z_{g,\min}} \\{1 - {\int_{z_{g,\max}}^{z}{{p\left( z_{g} \right)}{dz}_{g}}}} & {{{if}\mspace{14mu} z_{g,\max}} \leq z \leq z_{g,\min}} \\1 & {{{if}\mspace{14mu} z} \leq z_{g,\max}}\end{matrix}\quad} \right.$ wherein, the determination of the step /d/is a function of the value determined in the step /c-3/for each mesh cof dimension z of said column and for each facies A among the pluralityof facies, multiplied by the value determined in the step /c-4/ for eachmesh c of dimension z of said column.
 6. The method for determining amap of variance of height of liquid hydrocarbon in a geological model,the method comprises, for at least one column of the geological model:/e/ determining a value of expectation of height of liquid hydrocarbonfor said column, according to claim 1; /f/ determining a value ofvariance of height of liquid hydrocarbon for said column analyticallyaccording to the sum determined in the step /e/.
 7. The method accordingto claim 6, wherein, the values of the proportion cubes for twodifferent cells of the model and for different architectural elementsare considered independently between them.
 8. The method according toclaim 6, wherein, the method further comprises: /g/ receiving acorrelogram ρ(A,Δz) according to a direction of said column for saidfacies A in a given architectural element; /h/ determining theexpectation of the product of the presence of the facies A in saidarchitectural element in a mesh c1 by the presence of the facies A insaid architectural element in a mesh c2, with the distance between c1and c2 according to said direction of the correlogram being Δz, with theprobability of the presence of the facies A of said architecturalelement, for the mesh c1 being p_(A,AE,c1), with the probability of thepresence of the facies A of said architectural element for the mesh c2being p_(A,AE,c2), according to:(1−ρ(A,Δz))·(Esp(p _(A,AE,c1))·Esp(p _(A,AE,c2))+√{square root over(Var(p _(A,AE,c1))·Var(p _(A,AE,c2)))})+ρ(A,Δz)·(½Esp(p_(A,AE,c1))+½·Esp(p _(A,AE,c2))) wherein, the determination of the step/f/ is according to the determination of the step /h/.
 9. A device fordetermining a map of expectation of height of liquid hydrocarbon in ageological model, with the geological model comprising meshes able to beassociated with a facies among a plurality of facies and with anarchitectural element among a plurality of architectural elements, thedevice comprises, for at least one column of the geological model: /a/an interface for receiving an appointment cube of architectural elementfor each architectural element AE among the plurality of architecturalelements, with the apportionment cube describing for each mesh c of themodel the probability p_(AE)(c) that said mesh belongs to saidarchitectural element, /b/ an interface for the receiving a proportioncube for each facies among the plurality of facies and for eacharchitectural element among the plurality of architectural elements,with the proportion cube associated with a facies A and with anarchitectural element AE describing for each mesh c of the model theproportion of said facies in said mesh if said mesh belongs to saidarchitectural element AE, with said proportion being a random variableof triangular distribution defined by a minimum value p_(A,AE,min)(c), amaximum value p_(A,AE,max)(c), and a mode p_(A,AE,mode)(c), /c/ acircuit for the determining, for each mesh c of said column and for eachfacies A of the plurality of facies, the sum over the plurality ofarchitectural elements ofp_(AE)(c)·⅓(p_(A,AE,min)(c)+p_(A,AE,max)(c)+p_(A,AE,mode)(c)). /d/ acircuit for the determining of a value of expectation of height ofliquid hydrocarbon for said column analytically according to the sumdetermined by the circuit/c/ for each mesh c of said column and for eachfacies A.
 10. A computer program product comprising instructions for theimplementing of the method according to claim 1, when this program isexecuted by a processor.